{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "# General E(2)-Equivariant Steerable CNNs  -  Hands-on tutorial"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We start by importing the necessary packages. The user typically only needs to interact with the high level functionalities provided in the subpackages `e2cnn.gspaces` and `e2cnn.nn`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import torch\n",
    "\n",
    "from e2cnn import gspaces\n",
    "from e2cnn import nn"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Basic Examples\n",
    "\n",
    "Throughout the whole tutorial we consider data which is given as a signal on the plane $\\mathbb{R}^2$, for instance images.\n",
    "We formalize such signals as **feature vector fields**, i.e. functions\n",
    "$$ f: \\mathbb{R}^2 \\to \\mathbb{R}^c $$\n",
    "\n",
    "which assign a $c$-dimensional feature vector $f(x)\\in\\mathbb{R}^c$ to each spatial position $x\\in\\mathbb{R}^2$.\n",
    "\n",
    "The plane $\\mathbb{R}^2$ has many symmetries which can be exploited by equivariant CNNs.\n",
    "Our library focuses on exploiting **isometries**, i.e. distance preserving transformations.\n",
    "Specifically for $\\mathbb{R}^2$, the isometries are formalized by the **Euclidean group** $E(2)$ which consists of *translations*, *rotations* and *reflections*.\n",
    "\n",
    "E(2) steerable CNNs can adaptively choose the level of symmetries which they respect.\n",
    "Since we are choosing a *convolutional* network design, all models will be equivariant under translations.\n",
    "The open choice is therefore in the level of point symmetries (here reflections and rotations) which are being considered.\n",
    "All of these choices are subgroups $G\\leq O(2)$ of the orthogonal group.\n",
    "For simplicity, we will in the following consider a the cyclic subgroup $C_4$, which models the $4$ *rotations* which are multiples of $\\frac{\\pi}{2}$, that is, $\\big\\{0, \\frac{\\pi}{2}, \\pi, \\frac{3\\pi}{2}\\big\\}$.\n",
    "Because these are perfect symmetries of the grid, transforming an image with this group does not require any interpolation.\n",
    "\n",
    "We determine both the **point group** and its **action on the space** by instantiating a subclass of `gspace.GSpace`.\n",
    "For the rotational action of $G=C_4$ on $\\mathbb{R}^2$ this is done by:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "r2_act = gspaces.Rot2dOnR2(N=4)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Having specified the symmetry transformation on the *base space* $\\mathbb{R}^2$, we next need define how the signals sampled on this space transform under the action of $G$.\n",
    "This transformation law of feature fields is implemented as a **FieldType**.\n",
    "\n",
    "Each single feature space (layer) of E(2) steerable CNNs has its own transformation law which the user needs to specify.\n",
    "The transformation law of the network *input* and *output* is typically determined by the inference task.\n",
    "\n",
    "As a trivial, yet important example we consider the case of gray-scale images as network input.\n",
    "A rotation of a gray-scale image is performed by moving each pixel to a new position without changing their intensity values.\n",
    "The invariance of the scalar pixel values under rotations is modeled by the **trivial representation** $\\rho:G\\to\\operatorname{GL}(1),\\ g\\mapsto(1)$ of $G$ and identifies them as **scalar fields**.\n",
    "Formally, a scalar field is a function $f: \\mathbb{R}^2 \\to \\mathbb{R}$ mapping to a feature vector with $c=1$ channels.\n",
    "A rotation by an angle $\\theta \\in C_4$ transforms this scalar field as\n",
    "\n",
    "$$ \\big[\\mathcal{R}_{\\theta}\\, f\\big](x)\n",
    "   \\ :=\\ \\rho(\\theta)\\,f\\big(\\psi(-\\theta)x\\big)\n",
    "   \\ =\\ 1\\cdot f\\big(\\psi(-\\theta)x\\big)\n",
    "   \\ =\\ f\\big(\\psi(-\\theta)x\\big),$$\n",
    "\n",
    "where $\\mathcal{R}$ is the rotation operator acting on the field and $\\psi(\\theta)\\in SO(2)$ is a rotation matrix.\n",
    "\n",
    "We instantiate the `nn.FieldType` modeling a gray-scale image by passing it the `gspaces.GSpace` instance and the trivial representation.\n",
    "The latter is passed as a *list of representations* which allows for feature spaces that comprise multiple independent feature fields as used in the hidden layers later on."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "feat_type_in = nn.FieldType(r2_act, [r2_act.trivial_repr])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "When we build a model **equivariant** to a group $G$, we require that the output produced by the model transforms consistently when the input transforms under the action of an element $g \\in G$.\n",
    "For a function $F$ (e.g. a neural network), the **equivariance constraint** requires:\n",
    "\n",
    "$$ \\mathcal{T}^\\text{out}_g \\big[F(x)\\big]\\ =\\ F\\big(\\mathcal{T}^\\text{in}_g[x]\\big) \\quad \\forall g\\in G$$\n",
    "\n",
    "where $\\mathcal{T}^\\text{in}_g$ is the transformation of the input by the group element $g$ while $\\mathcal{T}^\\text{out}_g$ is the transformation of the output by the same element.\n",
    "We have already defined $\\mathcal{T}^\\text{in}$ by the *field type* `feat_type_in`. \n",
    "The transformation law $\\mathcal{T}^\\text{out}$ of the output of the first layer is similarly chosen by defining an instance `feat_type_out` of `nn.FieldType`.\n",
    "\n",
    "Instead of scalar feature fields, we often choose so called **regular feature fields** $f:\\mathbb{R}^2\\to\\mathbb{R}^{|G|}$ in the hidden layers of the model.\n",
    "These fields correspond to *group convolutions* and are empirically found to work best in most cases.\n",
    "Regular feature fields associate a feature vector of dimensionality equal to the cardinality of $G$, in our case $|C_4|=4$, to each point.\n",
    "Their name comes from the fact that they transform under the **regular representation** $\\rho_\\text{reg}: G \\to \\operatorname{GL}\\big(\\mathbb{R}^{|G|}\\big)$ of $G$ which acts via permutations.\n",
    "Formally, a regular feature field thus transforms according to:\n",
    "\n",
    "$$ \\big[\\mathcal{R}_{\\theta}\\, f\\big](x)\\ :=\\ \\rho_\\text{reg}(\\theta)\\, f\\big(\\psi(-\\theta)x\\big)$$\n",
    "\n",
    "A feature space consisting of a single regular feature field is instantiated as before, with the only difference that the trivial representation is replaced by the regular representation:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "feat_type_out = nn.FieldType(r2_act, [r2_act.regular_repr])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As deep feature spaces typically comprise multiple independent features (equal to channels for conventional CNNs), we usually choose multiple feature fields.\n",
    "For example, we can use 3 regular fields.\n",
    "The full feature space is in this case modeled as a *stacked* field $f: \\mathbb{R}^2 \\to \\mathbb{R}^{3|G|}$ which transforms according to the **direct sum** of three regular representations:\n",
    "\n",
    "$$\\rho(\\theta)\n",
    "    \\ =\\ \\rho_\\text{reg}(\\theta) \\oplus \\rho_\\text{reg}(\\theta) \\oplus \\rho_\\text{reg}(\\theta)\n",
    "    \\ =\\ \\begin{bmatrix} \n",
    "            \\rho_\\text{reg}(\\theta) & 0 & 0 \\\\\n",
    "            0 & \\rho_\\text{reg}(\\theta) & 0 \\\\\n",
    "            0 & 0 & \\rho_\\text{reg}(\\theta) \\\\\n",
    "          \\end{bmatrix}\n",
    "          \\quad\\in\\ \\mathbb{R}^{3N \\times 3N}$$\n",
    "\n",
    "Intuitively, the direct sum simply builds a block diagonal representation of the three regular representations, which implies that the three regular fields transform independently from each other.\n",
    "\n",
    "We instantiate a `nn.FieldType` composed of 3 regular representations by passing the full field representation as a list of three regular representations:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "feat_type_out = nn.FieldType(r2_act, 3*[r2_act.regular_repr])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Once having defined how the input and output feature spaces should transform, we can build neural network functions as **equivariant modules**.\n",
    "These are implemented as subclasses of an abstract base class `nn.EquivariantModule` which itself inherits from `torch.nn.Module`.\n",
    "\n",
    "We start by instantiating a convolutional layer that maps between fields of types `feat_type_in` and `feat_type_out`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "conv = nn.R2Conv(feat_type_in, feat_type_out, kernel_size=3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Each equivariant module has an input and output type.\n",
    "As a function (`.forward()`), it *requires* its inputs to transform according to its input type and is guaranteed to return feature fields associated with its output type.\n",
    "To prevent the user from accidentally feeding an incorrectly transforming input field into an equivariant module, we perform a dynamic type checking.\n",
    "In order to do so, we define **geometric tensors** as data containers.\n",
    "They are wrapping a *PyTorch* `torch.Tensor` to augment them with an instance of `FieldType`.\n",
    "\n",
    "Let's build a few random 32x32 gray-scale images and wrap them into an `nn.GeometricTensor`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "x = torch.randn(4, 1, 32, 32)\n",
    "x = nn.GeometricTensor(x, feat_type_in)\n",
    "\n",
    "assert isinstance(x.tensor, torch.Tensor)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As usually done in *PyTorch*, an image or feature map is stored in a 4-dimensional array of shape BxCxHxW, where B is the batch-size, C is the number of channels and W and H are the spatial dimensions."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can feed a geometric tensor to an equivariant module as we feed normal tensors in *PyTorch*'s modules:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "y = conv(x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can verify that the output is indeed associated with the output type of the convolutional layer:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "assert y.type == feat_type_out"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lets check whether the output transforms as described by the output type when the input transforms according to the input type.\n",
    "The $G$-transformation of a geometric tensor is hereby conveniently done by calling `nn.GeometricTensor.transform()`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# for each group element\n",
    "for g in r2_act.testing_elements:\n",
    "    # transform the input with the current group element according to the input type\n",
    "    x_transformed = x.transform(g)\n",
    "    \n",
    "    # feed the transformed input in the convolutional layer\n",
    "    y_from_x_transformed = conv(x_transformed)\n",
    "    \n",
    "    # the result should be equivalent to rotating the output produced in the \n",
    "    # previous block according to the output type\n",
    "    y_transformed_from_x = y.transform(g)\n",
    "    assert torch.allclose(y_from_x_transformed.tensor, y_transformed_from_x.tensor, atol=1e-5), g"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Any network operation is required to be equivariant.\n",
    "`e2cnn.nn` provides a wide range of equivariant network modules which guarantee this behavior.\n",
    "\n",
    "As an example, we will next apply an **equivariant nonlinearity** to the output feature field of the convolution.\n",
    "For the specific case of regular representations, any pointwise nonlinearity like *ReLUs* are equivariant.\n",
    "Note that this is *not* the case for many other choices of representations / field types!\n",
    "\n",
    "We instantiate a `e2cnn.nn.ReLU`, which, as an `nn.EquivariantModule`, requires to be informed about its input type to be able to perform the type checking.\n",
    "Here we are passing `feat_type_out`, the output of the equivariant convolution layer, as input type.\n",
    "It is not necessary to pass an output type to the nonlinearity since this is here determined by its input type."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "relu = nn.ReLU(feat_type_out)\n",
    "\n",
    "z = relu(y)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can verify the equivariance again:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "# for each group element\n",
    "for g in r2_act.testing_elements:\n",
    "    y_transformed = y.transform(g)\n",
    "    z_from_y_transformed = relu(y_transformed)\n",
    "    \n",
    "    z_transformed_from_y = z.transform(g)\n",
    "    \n",
    "    assert torch.allclose(z_from_y_transformed.tensor, z_transformed_from_y.tensor, atol=1e-5), g"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In *deep learning* we usually want to stack multiple layers to build a deep model.\n",
    "As long as each layer is equivariant and consecutive layers are compatible, the equivariance property is preserved by induction.\n",
    "\n",
    "The compatibility of two consecutive layers requires the output type of the first layer to be equal to the input type of the second layer.\n",
    "\n",
    "In case we feed an input with the wrong type to a module, an error is raised:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Error! the type of the input does not match the input type of this module\n"
     ]
    }
   ],
   "source": [
    "layer1 = nn.R2Conv(feat_type_in, feat_type_out, kernel_size=3)\n",
    "layer2 = nn.ReLU(feat_type_in) # the input type of the ReLU should be the output type of the convolution\n",
    "\n",
    "x = nn.GeometricTensor(torch.randn(3, 1, 7, 7), feat_type_in)\n",
    "\n",
    "try:\n",
    "    y = layer2(layer1(x))\n",
    "except AssertionError as e:\n",
    "    print(e)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Simple deeper architectures can be built using a **SequentialModule**:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "feat_type_in = nn.FieldType(r2_act, [r2_act.trivial_repr])\n",
    "feat_type_hid = nn.FieldType(r2_act, 8*[r2_act.regular_repr])\n",
    "feat_type_out = nn.FieldType(r2_act, 2*[r2_act.regular_repr])\n",
    "\n",
    "model = nn.SequentialModule(\n",
    "    nn.R2Conv(feat_type_in, feat_type_hid, kernel_size=3),\n",
    "    nn.InnerBatchNorm(feat_type_hid),\n",
    "    nn.ReLU(feat_type_hid),\n",
    "    nn.R2Conv(feat_type_hid, feat_type_hid, kernel_size=3),\n",
    "    nn.InnerBatchNorm(feat_type_hid),\n",
    "    nn.ReLU(feat_type_hid),\n",
    "    nn.R2Conv(feat_type_hid, feat_type_out, kernel_size=3),\n",
    ").eval()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As every layer is equivariant and consecutive layers are compatible, the whole model is equivariant."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "x = torch.randn(1, 1, 17, 17)\n",
    "x = nn.GeometricTensor(x, feat_type_in)\n",
    "\n",
    "y = model(x)\n",
    "\n",
    "# for each group element\n",
    "for g in r2_act.testing_elements:\n",
    "    x_transformed = x.transform(g)\n",
    "    y_from_x_transformed = model(x_transformed)\n",
    "    \n",
    "    y_transformed_from_x = y.transform(g)\n",
    "    \n",
    "    assert torch.allclose(y_from_x_transformed.tensor, y_transformed_from_x.tensor, atol=1e-5), g"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Usually, at the end of the model we want to produce a single feature vector to use for classification.\n",
    "To do so, it is common to pool over the spatial dimensions, e.g. through average pooling.\n",
    "\n",
    "This produces (approximatively) translation-invariant feature vectors."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "torch.Size([1, 8, 1, 1])\n"
     ]
    }
   ],
   "source": [
    "avgpool = nn.PointwiseAvgPool(feat_type_out, 11)\n",
    "\n",
    "y = avgpool(model(x))\n",
    "\n",
    "print(y.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In our case, the feature vectors $f(x)\\in\\mathbb{R}^c$ associated to each point $x\\in\\mathbb{R}^2$ have a well defined transformation law.\n",
    "The output of the model now transforms according to `feat_type_out` (here two $C_4$ regular fields, i.e. 8 channels).\n",
    "For our choice of regular representations (which are permutation representations) the channels in the feature vectors associated to each point permute when the input is rotated."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "rotation by 0*pi/2: [0.3197 0.3013 0.3031 0.3098 0.2218 0.2497 0.2616 0.2284]\n",
      "rotation by 1*pi/2: [0.3098 0.3197 0.3013 0.3031 0.2284 0.2218 0.2497 0.2616]\n",
      "rotation by 2*pi/2: [0.3031 0.3098 0.3197 0.3013 0.2616 0.2284 0.2218 0.2497]\n",
      "rotation by 3*pi/2: [0.3013 0.3031 0.3098 0.3197 0.2497 0.2616 0.2284 0.2218]\n"
     ]
    }
   ],
   "source": [
    "for i in range(4):\n",
    "    print(f'rotation by {i}*pi/2:', y.transform(i).tensor[0, ...].detach().numpy().squeeze())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Many learning tasks require to build models which are **invariant** under rotations.\n",
    "We can compute invariant features from the output of the model using an **invariant map**.\n",
    "For instance, we can take the maximum value within each regular field.\n",
    "We do so using `nn.GroupPooling`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "rotation by 0*pi/2: [0.3197 0.2616]\n",
      "rotation by 1*pi/2: [0.3197 0.2616]\n",
      "rotation by 2*pi/2: [0.3197 0.2616]\n",
      "rotation by 3*pi/2: [0.3197 0.2616]\n"
     ]
    }
   ],
   "source": [
    "invariant_map = nn.GroupPooling(feat_type_out)\n",
    "\n",
    "y = invariant_map(avgpool(model(x)))\n",
    "\n",
    "for i in range(4):\n",
    "    print(f'rotation by {i}*pi/2:', y.transform(i).tensor[0, ...].detach().numpy().squeeze())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "# for each group element\n",
    "for g in r2_act.testing_elements:\n",
    "    # rotated the input image\n",
    "    x_transformed = x.transform(g)\n",
    "    y_from_x_transformed = invariant_map(avgpool(model(x_transformed)))\n",
    "    \n",
    "    y_transformed_from_x = y # no .transform(g) needed since y should be invariant!\n",
    "    \n",
    "    # check that the output did not change\n",
    "    # note that here we are not rotating the original output y as before\n",
    "    assert torch.allclose(y_from_x_transformed.tensor, y_transformed_from_x.tensor, atol=1e-6), g"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Other Field Field Types"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Besides scalar fields and regular fields, many other types of feature fields exist.\n",
    "\n",
    "An example of practical importance are (tangent) **vector fields** which model for instance an optical flow fields.\n",
    "On the plane, a vector field\n",
    "$$f: \\mathbb{R}^2 \\to \\mathbb{R}^2$$\n",
    "\n",
    "associates a 2-dimensional (tangent) vector to each position.\n",
    "Under rotations a vector field transforms by:\n",
    "1. moving each vector to a new location\n",
    "2. rotating the vectors themself\n",
    "\n",
    "Together, the $G$-transformation of a vector field is therefore given by\n",
    "\n",
    "$$ \\big[\\mathcal{R}_{\\theta}\\, f\\big](x)\n",
    "   \\ :=\\ \\psi(\\theta\\,)\\, f\\big(\\psi(-\\theta)x\\big),$$\n",
    "\n",
    "that is, the representation $\\rho(\\theta)$ coincides with the rotation matrix\n",
    "$\\psi(\\theta) = \\begin{bmatrix} \\cos(\\theta) & \\sin(\\theta) \\\\ -\\sin(\\theta) & \\cos(\\theta) \\end{bmatrix} \\in SO(2)$ itself.\n",
    "\n",
    "In general, the *type* of a feature field is in one-to-one correspondence with a choice of group representation $\\rho:G\\to\\operatorname{GL}(\\mathbb{R}^c)$ which describes how the feature vectors transform under the action of $G$.\n",
    "`e2cnn` implements **irreducible representations** (irreps), **regular representations**, **quotient representations** and **induced representations** and provides an interface to construct any other choice of representation.\n",
    "For more information we refer to the\n",
    "[docs](https://quva-lab.github.io/e2cnn/api/e2cnn.group.html#representations)\n",
    "and to Section 2.6 of our\n",
    "[paper](https://arxiv.org/pdf/1911.08251.pdf)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For now, let us return to the example of (tangent) vector fields.\n",
    "The rotation matrix representation $\\psi$ of the cyclic group $G=C_4$ is identical to the *irreducible representation* (irrep) of $C_4$ with frequency $1$.\n",
    "To specify a vector field as output, we therefore have to pass this irrep to the `FieldType`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "feat_type_out = nn.FieldType(r2_act, [r2_act.irrep(1)])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To verify that the output field indeed transforms as a *vector field*, we build a simple model which maps a random scalar field to a vector field and apply it to rotated versions of the scalar input field.\n",
    "As expected, the output vector field transforms by moving each vector to a new position and rotating them."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA6IAAAEICAYAAABF4EkmAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOzdd3hU1fo24GelUwISeg8KVgxikCLogQihBalSpQnCD4PoAUEUCeVI1XNEBalGOhICAenSCUSMSSiRJjWUAClACoGUmff7IyEfJSSTkMyaSZ77uuZyZu81ex7GrJl599p7bSUiICIiIiIiIjIXG90BiIiIiIiIqGhhIUpERERERERmxUKUiIiIiIiIzIqFKBEREREREZkVC1EiIiIiIiIyKxaiREREREREZFYsRImIiIiIiMisWIgWQUopF6VUgFLqjlIqQinVW3cmoqJMKfWSUmq3UipOKXVWKdX5gXXvKKVOKaWSlFJ7lFI1dWYlKgqUUsOVUiFKqWSl1OJH1nVXSp1USiUopU4opTo9sv7fSqnrGf3ZVynlaNbwRIWcUspRKfVzxm/YBKXUYaVU24x1fZRSiQ/ckpRSopRyz1ivlFIzlFKxGbeZSiml919UdLEQLZrmAEgBUBFAHwBzlVKv6I1EVDQppewAbACwCYALgCEAliulnldKlQOwDsD4jHUhAFbrykpUhEQC+BqA74MLlVJVASwHMBJAKQCjAaxUSlXIWN8awFgA7wBwBfAsgElmS01UNNgBuAzgXwBKI/070k8p5SoiK0Sk5P0bgI8AnAcQlvHcIQA6AagHwA2AF4Ch5v4HUDolIrozkBkppUoAuAWgroj8k7FsGYCrIjJWaziiIkgpVRfAIQDOkvGBrJT6HcCfSP+iHSAib2YsLwEgBkB9ETmlKTJRkaGU+hpANREZkPG4EYCNIlLhgTbRAN4VkT+UUisBXBSRLzPWvQNghYhUMn96oqJDKXUMwCQRWfvI8j0A9orIpIzHQQAWi8iCjMeDAHwoIo3NnZk4IloUPQ/AcL8IzXAUAEdEifTI6pAgBaAu0vvl0fsLReQOgHNgfyXSJQTASaXUu0op24zDcpMBHMtY/1CfzbhfUSlV1sw5iYoMpVRFpP++Pf7I8poA3gaw9IHFWfVRfqdqwkK06CkJIO6RZXEAnDVkISLgFIAoAKOVUvZKKU+kH25UHOyvRBZFRAxI/1G7EukF6EoAQzN2EgGP99n799lniQqAUsoewAoAS7I4UqgfgEARufDAsqz6aEmeJ6oHC9GiJxHp57U8qBSABA1ZiIo8EUlF+vkq7QFcBzAKgB+AK2B/JbIoSqmWAGYCaA7AAek7jRYppV7LaPJon71/n32WKJ8ppWwALEP6vCfDs2jSD8CSR5Zl1UcThecqasFCtOj5B4CdUqrOA8vq4ZHDGYjIfETkmIj8S0TKikhrpE9wEoz0flnvfruMc0SfA/srkS6vAdgvIiEiYhSRv5B+PnfLjPUP9dmM+zdEJNbMOYkKtYwRzJ+RPvFm14ydug+ubwqgCgD/R56aVR/ld6omLESLmIzDh9YBmKyUKpHRUTsifY8SEWmglHJTSjkppYorpT4DUBnAYgABAOoqpboqpZwA+AA4xomKiAqWUsouo8/ZArDN6J92AP4C8Nb9EVClVH0Ab+H/nyO6FMAgpdTLSqkyAL5Cel8movw1F8BLADqIyN0s1vcHsFZEHj0aYSmAkUqpqkqpKkg/CmlxgSalJ2IhWjR9BKAY0s9LWwVgmIhwbxCRPn0BXEN6n3wHQCsRSRaRaABdAUxB+mzXjQD01JaSqOj4CsBdpF+K5f2M+1+JyD4AEwH4K6USAKwFMFVEfgcAEdmG9EN39wCIyLhNMHt6okIsYxKioUg/QuH6A9cM7ZOx3glAdzx+WC4AzAewEUA4gL8BbM5YRhrw8i1ERERERERkVhwRJSIiIiIiIrNiIUpERERERERmxUKUiIiIiIiIzIqFKBEREREREZmVna4XLleunLi6uup6eSKLEBoaGiMi5XXnyAr7KBH7KJnm3LlzuHv3LurWras7SpHDPkpk2bLro9oKUVdXV4SEhOh6eSKLoJSK0J3hSdhHidhHKWdJSUkoV64cPv30U0yfPl13nCKHfZTIsmXXR3loLhEREVEe7dixA3fv3kXHjh11R8lRdHS07ghERJlYiBIRERHl0YYNG1CxYkU0atRId5QcTZ48GXfu3NEdg4gIAAtRIiIiojwxGAzYuHEjOnbsCBsby/5JdefOHSxZsgSHDx/WHYXI7K5cuaI7wmNWrFiBS5cu6Y6hlWV/apJZrF27VncEKqTS0tJ0RyAiKjBBQUGIiYmxisNy/fz8kJCQgODgYN1RiMzq77//xrhx43THeEzt2rXh5uaGX3/9VXcUbViIFnHff/89Vq5cqTsGWbn4+Hj8+eef+OWXXzB69Gh4eXmhfv36CAoK0h2NiKjArF+/HiVLloSHh4fuKDlatGgRALAQpSLlzJkzaNmyJapWrao7ymMaNmyI6tWro1evXujbty/i4uJ0RzI7k2fNVUrZAggBcFVEvB5Z5whgKQB3ALEAeojIxXzMCQAQERw9ehQvvfQSHB0d83vzRc5PP/2ETz/9FF9//bXuKGTFDh06hPbt2+PmzZuZy8qXL49t27bh9ddf15iMiKjgiAjWr1+Ptm3bwsnJSXecbB0/fjxzxyALUcuSkJCAuXPn4vbt24iPj3/o1rBhQ0ycONHi/74sVUREBN555x3cuHED9erV0x3nMUopDB48GJ9++imWL1+OwMBALF++HM2aNdMdzWxyMyL6CYCTT1g3CMAtEakN4DsAM542WFaUUvD390epUqXQpEkT/Pvf/8bq1atx8eJFiEhBvGShtWDBAnh7ewMAiwXKM4PBgMuXL6N8+f9/eaiaNWviwIED/Lt6SgaDAb169ULNmjXx0ksvwd3dHW+99RY8PT3x8ccf49atW7ojkoW7efMm5s2bhyNHjvAw+QJw/PhxnD9/3ioOy70/GgoAFy5c4Oy5FsTZ2Rnx8fGYNm0a5syZg2XLlmHjxo1o2bIlpk6dyiI0j65du4Z33nkHly9fBgC4ublpTpS1999/Hw4ODgCAq1evYvHixUWqf5o0IqqUqgagPYApAEZm0aQjgIkZ9/0BzFZKKXmK6tBgMGD8+PGIiYlBbGwsYmJiEBMTg6ioKKSkpODQoUM4dOgQAMDR0RHe3t6YMGECSpUqldeXfGp3795FsWLFkJaWBjs7bZdozVFCQsJDJ0dbS8Fw9+5dTJs2DZMmTYJSSnecHIkIRMTiJ7DIq5MnT6JLly44deoU6tevjxo1asDZ2Rnbt2+3yENgAODy5cuoXr267hgmsbGxgY2NzUN91c7ODiNHjsT48eNRsmRJjemydvbsWdSsWRP29va6oxAAo9GImTNn4sKFCyhevDjeeOMNNG7cGB06dEDTpk11x8vS1atXLfbz41Hr16+HnZ0d2rVrpztKjlq2bAlbW1vMnTsXP/zwA86dO/fQDkRLJSJW8X2fV/v378cnn3yCI0eOZC6rVq0aVq1aZbGjYufPn8ezzz6rO0aOLl68iHr16uHcuXNwdHREnTp1dEfKUtmyZdGlSxds374dt27dwttvv20VfTPf3P+xnN0N6cWlO4DmADZlsf5vANUeeHwOQLks2g1B+uG9ITVq1JCcVK5cWZ577jlp1KiRtGvXTvr16yf9+vUTAAJAXn/9dZk9e7bcvHkzx20VtNDQUKlYsaIEBgaKt7e37jg5Sk1NlfXr18v48eN1RzFZt27dBIBMmDBBd5QcJSQkSOPGjWX+/PnZtgMQIib0QXPdctNH79y5I15eXrJ582YxGo0yevRoiY2NNfk9Mrfjx4+Lk5OTjBkzRlJTU3XHMUmHDh3kmWeeEQDSvHlzOX78uO5IT5SSkiLPPvus1K1bVw4ePKg7Trb+/vtvWbx4sTRt2lQmT54s0dHRT2xrzX303Llz4uLiIgDExsZG2rVrJ/7+/pKcnJyr98tc9u7dKw4ODjJnzhzdUUwSFRUlmzdv1h3DZBEREbJv3z7dMUwWGhoqL774ovzzzz/ZtrPmPnrs2DFp0aKF+Pr6ioeHh3h5eUlMTIzJ75G5BQQEiJ2dnUyZMkXS0tJ0x8lRbGys/Pbbb/LJJ5/ojpKtnTt3ypYtW+S7776T27dv645jMl9fX2nVqlWO7bLro6Z0KC8AP2Xcf1IhejyLQrRsdtt1d3fP0z96wYIFMnz4cDl8+HCenl9Q7ty5Iy+++KIUK1ZMbG1ts/1hY0mMRqPuCCaJj4+XQYMGiaurqwCQyZMn646UIw8PD3F1dZWUlBS5detWlm0s7Qv0wVtu+6il/y2lpqbK2LFjBYC0aNFCbty4IXv37pUbN27ojpatRo0ayYoVKyz+/RUR2bVrlzz//PMCQD788EOL3TFx8+ZNcXR0FABSrVq1bNtacx81Go3SvXt3mTJlily+fNm0N0ejpKQk6dmzpwCQYcOGSUpKiu5IpFFMTIzY29vLhAkTxN/f/4ntrLmPPmjLli0W/zmfmJgoQ4YMEQDy1ltvycWLF3VHMonBYNAdIVsGg8Hi/99nZebMmaKUynGnxNMWotMAXAFwEcB1AEkAlj/SZjuAJhn37QDEAFDZbTevhagl/o8yGo3yxRdfSOnSpTNHaxctWqQ7VqF1/fp12bhxo1y7dk13lGzt3r1bAMi4ceOkZ8+eWbYpLF+g1mTdunXi7OwsVatWFS8vL/Hy8rLIzxWR9C8na9o7KiJy7949mTRpkjg6Okr58uVl2bJlFvn+9u7dWwBIt27dsm1nzX3UaDRa5HufHaPRKFOmTMncYWTJo0NUcFJSUmT//v3y1ltviYuLi9ja2j7xSBZr7qPWKiAgQMqWLSulS5eWlStXiohIWFiYxRd8lL9+/vlnAZDj4Ft2fTTHk9dE5AsRqSYirgB6AtgtIu8/0uw3AP0z7nfLaFMgswdZ4rkCSilMnDgR/fr1y1y2Zs0ajYkKt4oVK8LLywuVKlXSHeWJbt68icOHD6NMmTKYMmVK5vnMpF/nzp0REhKCEiVKYNOmTdi0aRPmz5+vO1aWbGxsULp0ad0xcsXR0RE+Pj4IDw+Hm5sb+vbti1atWuGff/7RHe0hgwYNAgA0btxYc5KCo5SyyO/M7Cil8OWXXyIgIADBwcFo1KgRTpw4oTsWmZm9vT38/f0RGBiImzdvwmAwICEhQXcsytCpUyccO3YMDRs2RO/evdG3b198++23mDZtmu5oZEZly5YFAMTExOR5G3meRUUpNVkp9W7Gw58BlFVKnUX6ZEZj85zISjk4OOCHH37A0qVL4eTkhF27dj10OQsqWlxcXJCUlJQ5s2lERATu3LmjORXdFxgYiAsXLmQ+HjlyJE6dOqUxUeFTp04d7NixA8uXL0d4eDheffVVTJo0CcnJyZltCmh/pUmaN2+OWrVqFepC1Jp16tQJQUFBSEtLQ+PGjbFlyxbdkcjMZs6c+dBkivHx8RrT0KOqVKmCbdu24X//+x/8/PywcuVKjB8/Hr///rvuaGQm5cqVAwDExsbmeRu5KkRFZK9kXENURHxE5LeM+/dE5D0RqS0iDUXkfJ4TWbm+ffsiKCgI1apVw4YNG3THIY3GjRuHESNGAEj/wX369GnNiei+QYMG4ejRoxg1ahTKlSuHu3fvok+fPkhJSdEdrVBRSqFPnz44deoUBgwYgIkTJ8LNzQ179uyBwWDAlClTtGWzsbHB0KFDrWbW8KLIzc0NwcHBqFevHry8vPDtt99q3XlB5uXo6IjVq1fD2dkZAAtRS2RjY4PnnnsOlStXBpD+W6d3796IiIjQnIwK2okTJzIvOxMeHo4VK1bkaTuF87oSmtWvXx8hISH80CzilFL47rvv0KdPHwDplzshy/HSSy/h22+/xdWrV+Hv748KFSpg0qRJumMVSmXKlMH8+fNx8OBBODg4wMPDA++99x7Gjx8PX19fbblGjBiBYsWKaXt9ylmFChWwa9cufPDBBxg9ejQGDhz40Kg6FW61a9fGggULALAQtVRNmjTBF198gRYtWsDGxgaxsbHo1q0b7t27pzsaFaCYmJjMSwwNGzYMV65cydN2WIgWkLJly2aOhlHRZWNjg19++QXt2rXjeU4WysHBAV27dsXWrVsxbNgwGAwG3ZEKrTfffBNhYWGYMWMG1q9fDwAYMmQItm3bpiUPi1Dr4ODggIULF2LWrFlYtmwZWrRogevXr+uORWbSs2dPDBkyhIWohSpfvjyGDh2K3bt34+rVq5g9ezaKFSuGf//737qjUQF66623UKNGjczHeT3NhYVoAbK2SSKoYNjb22PNmjVwcXHRHYVyUK1aNdja2uqOUagZDAaEhIRkHmJpMBjQrVs3hIWFaU5GlkwphU8++QRbtmzBiRMn0LBhQxw+fBgAcP58kT0bqMiYNWsWqlatqjsG5aBSpUrw9vbG/v378dVXX+Hu3bu6I1EBUUrhgw8+AJA+6NKgQYM8bYeFKJEZFC9enHsHiQA4OTlh9erVOHr0KL755ht4enrCYDCgXbt2uHjxou54ZOFat26NP//8E05OTmjWrBnWrFmDNm3aICoqSnc0KkDFihWDm5ub7hiUC1WrVuVRJ4Vc//79YWNjAzc3N5QoUSJP22AhSmQmNjbsbkRA+p5UNzc3fPbZZ9i+fTtu3bqFFStWYO/evZyMhnL0wgsv4M8//8Sbb76J7t2748yZM/D29tYdi4ioSKlSpQratWv3VLPP2+VjHiIiolxzcnLCO++8ozsGWZG7d++iSpUqmY/9/f3h5+eH7t27a0xFRFS0fPDBB091/jaHaIiIiMiqVK5cGR9++CF69OgBO7v0fere3t48RJeIyIy8vLzQsmXLPD+fhSgRERFZFaUUmjVrhl9//RURERGYMGEC7Ozs8NFHH/HwbiIiM7G3t3+qicRYiBIREZHVqlKlCiZOnIiIiAh069YNZ86c0R2JiIhMwHNEiYiIyOo5ODigZ8+eumMQEZGJOCJKREREREREZsVClIiIiIiIiMyKhSgRERERERGZFQtRIiIiIiIiMisWokRERERERGRWLESJiIiIiIjIrHIsRJVSTkqpYKXUUaXUcaXUpCzaDFBKRSuljmTcBhdMXCIiIiIiIrJ2plxHNBmAh4gkKqXsARxQSm0VkUOPtFstIsPzPyIREREREREVJjkWoiIiABIzHtpn3KQgQxEREREREVHhZdI5okopW6XUEQBRAHaIyJ9ZNOuqlDqmlPJXSlV/wnaGKKVClFIh0dHRTxGbiAoC+yiRZWMfJbJs7KNEpjOpEBURg4i8BqAagIZKqbqPNNkIwFVE3ADsBLDkCdtZICINRKRB+fLlnyY3ERUA9lEiy8Y+SmTZ2EeJTJerWXNF5DaAvQDaPLI8VkSSMx4uBOCeL+mIiIiIiIio0DFl1tzySqlnMu4XA9ASwKlH2lR+4OG7AE7mZ0giIiIiIiIqPEyZNbcygCVKKVukF65+IrJJKTUZQIiI/AZghFLqXQBpAG4CGFBQgYmIiIiIiMi6mTJr7jEA9bNY7vPA/S8AfJG/0YiIiIgsT0pKCtLS0lC8eHHdUYiIrFauzhElIiIiKupCQ0OxYcMG3TGIiKwaC1EiIqIizmAw6I5gVfbv34/ly5frjkFEZNVYiBIRERVRIoKpU6ciKipKdxSrEhgYiO3bt+PGjRu6oxARWS0WokREREVQWloahg0bhpUrV6Jy5co5P4EApI8eHzhwAAaDAatXr9Ydh4jIarEQJSIiKmLu3LmDLl26YP78+WjVqpXuOFYlPDwccXFxAMDDc4mIngILUSIioiIkKioKHh4e2LhxIwDA09NTcyLrsn///sz7f/31F06fPq0xDRGR9WIhSkREVERERUXB09MTwcHBAAAHBwe8/fbbmlNZl5SUFIwaNQoAsHLlSly7dk1zIiIi68RClIiIqIioUKEC1q5dC0dHRzg4OODNN99EiRIldMeyKp999hk8PDwAAK+++iqaN2+uN9Ajrl+/juPHj+uOQUSUIxaiRERERYSIwNvbGyVLlsTu3bvh5eWlO5JVKlWqFAAgPj5ec5LHVahQAR07dsSRI0d0RyEiypad7gBERERkHn5+fti+fTt8fX3RtGlTNGrUSHckq2TJhaiNjQ1q1aoFDw8PbN++HW+88YbuSEREWeKIKAEADh8+rDsCEREVoNu3b+PTTz/F22+/jQEDBgAA7Oy4Pzov7heiCQkJmpNkrUmTJrh16xZatmyJoKAg3XGIiLLEQpRw+vRpjBw5UncMIiIqQOPGjUNsbCzmzZsHpZTuOFbNkkdEAaBx48YA0vN5enpi3759mhMRET2OhWgRFx0djXbt2qF06dK6oxARUQEJDg7G3LlzMWbMGLz00ku641g9Z2dnAJZfiAJAuXLlEBgYiOTkZI2JiIgex2NyirB79+6hU6dOOH/+PN59913dcYiIqACkpaVh6NChePbZZzFu3DjdcQoFe3t7FCtWzGILURcXF9StWxeurq7Yvn07evbsCUdHR92xiPLs9u3b+Prrr2E0GlGpUiVUrFjxof9WqFABNjYcX7M2LESLKKPRiAEDBmSeO+Lq6qo3EFmtpKQkGAyGzBECIrIsP/zwA44cOYLt27ejWLFiuuMUGqVKlbLYQhRIv8Zp1apVUbt2bYwdOxb+/v66IxHl2TPPPIORI0eia9euOHTo0EPrRo0ahZkzZ2pKRk+Duw6KKB8fH/j5+WU+rlWrlsY0ZK0SExPRrVs3ODk56Y5i8UQE0dHRCA0NRUBAAL7//nuMHDkSv//+u+5oVIhdunQJPj4+6NmzJzw9PXXHKVQsvRB99dVX4eLiAh8fH6xduxYHDhzQHYnoqVSpUgV79+7F4MGDH1p+9OhRBAQEIDU1VVMyyqscR0SVUk4A9gNwzGjvLyITHmnjCGApAHcAsQB6iMjFfE9L+WbChAmIiYnBihUrYGdnx0KUci0hIQHt27eHwWCAvb297jgWLSgoCP/3f/+H8PDwh5aXKlUKbdu2RWxsLMqWLaspHRVmI0aMgJ2dHb777jvdUQodSy9E7/voo48we/ZsjBo1Cn/88QcPX9TMaDTizJkzCA4Oxttvv42aNWvqjmRVHB0dsXDhQjRo0AAff/wxmjdvjkuXLqFbt26oVKkSPvjgA3z44YcWfaTf9evXcfXqVaSkpCA5ORkpKSmZ919++eUidR6/KZ9GyQA8RKQegNcAtFFKNX6kzSAAt0SkNoDvAMzIa6C4uDgcP348r08nE8XFxWHJkiX48MMPsXr1aovusGR5EhIS0KZNGwQGBqJhw4a645hMRLS87ptvvonQ0FD4+vo+tNPn/oyW5cqVQ5UqVdC6dWt89tlnWLJkCcLCwpCSkqIlLxUOGzZswIYNGzBt2jRUqlRJdxyTPNhHLfXSKPdZSyHq4OCAGTNmIDg4+KEjocg8IiMjsX79eowbNw6tWrWCi4sLXnzxRUydOpWHyj+FoUOHYs+ePejZsydOnjyJvXv3okWLFvjmm2/w7LPPol27dtiwYQPS0tJ0R32Ms7Mzfv31V7z11lto0aIFWrdujQ4dOmDkyJGoXLmy7ng5io2Nxc2bN/NlWzkWopIuMeOhfcbt0V9zHQEsybjvD+Adlce54YcPH445c+bgvffes4oPeAD48ccfsXXrVgQHB+uOYrIbN27A3d0dI0aMgKenp8Wf3/fnn3+iRo0aMBqNuqOYZO3atQgICNAdo8A4Oztj5syZaNq0KRo1aqQ7jknGjx+Pjz76SNvr29vbY+DAgTh9+jQWLVoEV1dXdO/eHQcOHMDcuXPRuXNn3L17Fz///DMGDBgAd3d3nD59Wlve3Nq4caNVHfo3b9487NmzR3eMAhUdHQ0PDw8MHTpUdxSTJCQkwM3NDVOnTkVUVBTGjh2LNWvW6I71RC1btkSTJk10xzBJly5d0K5dO6v5XQUAP/30E0JDQ3XHeGqXL1/Gpk2b8P3332Pnzp2Ii4sDAJw6dQrdu3fXnM40qampmD17Nrp27apth25WmjZtioEDB0IphX/9619YuXIlrl69ihkzZuDMmTPo1asXEhMTc96QmZUoUQLffPMNQkNDH/oNdenSJYvfAWcwGNC7d29s3LgR7777Li5evPh0GxSRHG8AbAEcAZAIYEYW6/8GUO2Bx+cAlMui3RAAIQBCatSoIY9atWqVIL3IlUqVKsmxY8cea2OJnn/+ebGxsZHnnntORowYIampqbojFTrJyckSGBioO4ZJ/v77bylRooS0bt1ajEZjtm0BhIgJfdBct5z66KMMBoPEx8eb9L7odunSJTlz5ozuGJlSUlJk48aNjy03Go1y+fJl2bJliyQnJ2tIlntpaWlSs2ZN6dChg+4oJtmyZYvY2NjI4MGDc2xr7X00p88gSxIRESGdOnUSGxsbsbe3l0qVKgkA8fHxEYPBoDue1bOmv4UDBw6IjY2NjBgxIse21tJH4+LiZP78+dKgQQMBIK1bt5Y9e/bk9S0qUFFRUSKS/jezfv16ef755wWAdOzYURISEjSnM43BYJATJ07ojpGjtLQ0mTNnjpQqVUqmT5+uO06Ovvzyy8xazcXFRXbv3p3jc7Lro7ntXM8A2AOg7iPLj2dRiJbNblvu7u4PhYyIiJDSpUs/VIiuWrXqKd4q86lbt25m7tWrV+uOQxoYjUZJSEiQuLg4ef7556VmzZoSExOT4/Ms7Qv0wdujfZQoKyEhIfLbb78JANm2bZvuODk6ceKElCpVSt588025d+9eju3ZR80vIiJCBgwYkPm9CkA6d+5sNT+A6encvn1bXF1d5ZVXXpGkpKQc21tjHz18+LBMnjzZIncO7N+/X7p06SIhISHSvHlzASDu7u6yd+9e3dEKtcjISAkKCtIdI1vr1q176HO5V69emTstspNdH83V5VtE5LZSai+ANhmjoPddAVAdwBWllB2A0gBMPnjYYDCgf//+iIuLg62tLTp06IDBgwejdevWuYmnzf1rc7Vt2xbvvfee5jSkQ0hICDZt2oTw8HBERETgwIEDnHyGigQvLy+UL18etWvXRt6EaVsAACAASURBVKtWrXTHydbNmzfx7rvv4plnnsG6det4XUUL5ezsjMjISDzzzDO4ffs2ACAgIABNmzbFhg0bOKdBITd8+HBERkbir7/+KrTnUL722mt47bXXdMd4zKFDh9CuXTsYDAasW7cO1apVw7Jly9C7d29OclXAKleubNHnh546dQr9+vUDANja2sLDwwMtWrRAHs/EzGTKrLnlAaRmFKHFALTE45MR/QagP4A/AHQDsDujAjbJf//7X1y5cgXTp09H//79rWZShfscHBxQrFgxzJkz56n/h5B1WrRoERYtWgSj0YhFixahQYMGuiMRFbjExERcv34d169fR8mSJeHl5YWAgACLLPBSU1PRvXt3REZG4uDBg6hYsaLuSPQEZcqUwfbt2wEAd+7cQWRkJK5evYrIyEjs27cP1apVg50dL4NemBiNRtjY2GDlypVYvnw5vvvuO7i5uemOVaSEhoaiTZs2medUvvrqqzh48KDFzyFCBS8+Ph7du3dH06ZN8d5776FTp075Nthiyid5ZQBLlFK2SJ/cyE9ENimlJiN9qPU3AD8DWKaUOov0kdCepgZIS0vDW2+9hdGjR1ttEefo6IgJEybwEihFVGJiIlauXJk5kdKcOXPwr3/9C7Vr19acjKhgRUZGZt4vUaIEFixYYJFFKACMHDkSu3btwtq1ay1yJIKyVqJECdSpUwd16tTRHYUK0MSJE9GrVy8MGzYMnp6eGDFihO5IRcqxY8fg6emZOZGSi4sLqlevjoMHD6JNmzaa05FuN2/exN69e+Hi4pLv286xEBWRYwDqZ7Hc54H79wDk6ZhUOzs7q5l17knc3d0xcuRI3TFIkzVr1mTuQXz99dexYsUKFqFUJFy9ehVA+lEhAQEBqFatmuZEWZs3bx5mz56NyZMno0uXLrrjENEDUlJSMHv2bHz//fewt7fH4sWLeRioGf3zzz/47LPP0L17dzRp0gSNGzdGnTp1rHZwiPJfQZ4OwWNb8sHEiRNhb2+vOwZpsmjRIiilMGbMGEyePBkODg66IxGZxf1CdMGCBRa7Q3Hv3r34+OOP0aNHD3z11Ve64xDRI3bu3Ilbt24BSB+JmzlzJr755hsefm0mtWvXxu+//647BhVR7OX5oGTJkrojkCYnT57E5cuXsXv3bjRv3lx3HCKzunr1KkaOHIn+/fvrjpKl8+fPo2vXrqhXrx58fX25h5/IAq1evTrzftOmTTFhwgQWoWbE0WfSiT2d6CnExsbi6NGjKFOmjO4oRGb36quvwtPTU3eMLMXHx6NDhw5wcHDAhg0bULx4cd2RiOgRycnJWL9+PWxsbPCf//wHY8eOZWFEVISwECV6Cs2aNdMdgUibdu3a6Y6QJYPBgN69e+PcuXPYt28fqlatqjsSEWVh+/btsLe3x7Zt2yz+8k9ElP+424mIiAqFAwcOIDU1FV9++SU2b96Mn3/+GY0aNdIdi4ie4OLFiwgLC2MRSlREcUSUiIgKBR8fH5QtWxb+/v4YO3Ys+vTpozsSEWXj448/5rnbREUYC1EiIrJ6x48fx549ewAAtWrVwuDBgyEi/JFLZMHYP4mKNh6aS0REVm/OnDmZ969evQp/f38YjUaNiYiIiCg7HBElIiKrFhcXh6VLlwIAGjZsiF9++QUvv/yy5lRERESUHY6IEhGRVVuyZAlSU1Mxffp0HDx4kEUoERGRFeCIKBERWS2j0YiQkBAcPnyYBSgREZEVYSFKRERWy2AwwNfXF3Z2/DojIiKyJvzmJiIiq2Vvb687AhEREeUBzxElIiIiIiIis2IhSkRERERERGaVYyGqlKqulNqjlDqplDqulPokizbNlVJxSqkjGTefgolLRERERERE1s6Uc0TTAIwSkTCllDOAUKXUDhE58Ui7QBHxyv+IREREREREVJjkOCIqItdEJCzjfgKAkwCqFnQwIiIiIiIiKpxydY6oUsoVQH0Af2axuolS6qhSaqtS6pUnPH+IUipEKRUSHR2d67BEVLDYR4ksG/sokWVjHyUyncmFqFKqJIC1AD4VkfhHVocBqCki9QD8CGB9VtsQkQUi0kBEGpQvXz6vmYmogLCPElk29lEiy8Y+SmQ6kwpRpZQ90ovQFSKy7tH1IhIvIokZ97cAsFdKlcvXpERERERERFQomDJrrgLwM4CTIvK/J7SplNEOSqmGGduNzc+gREREREREVDiYMmtuUwB9AYQrpY5kLPsSQA0AEJF5ALoBGKaUSgNwF0BPEZECyEtERERERERWLsdCVEQOAFA5tJkNYHZ+hSIiIiIiIqLCK1ez5hIRERERERE9LRaiREREREREZFYsRImIiIiIiMisWIgSERERERGRWbEQJSIiIiIiIrNiIUpWZ+fOnbojEBERFQr37t3THYGIiigWomRVUlJSMGzYMBiNRt1RiIiIrNqZM2fg6+urOwYRFVEsRMmq/PLLLzh79iyuX7+uOwoREZHVSkhIQKdOneDo6Kg7ChEVUSxEyWokJydjypQpAIBLly5pTkNERGSdRAQDBgzAiRMnUKlSJd1xiKiIYiFKVuOXX37B5cuXAbAQJSIiyqupU6di3bp1AICKFStqTkNERRULUbIKD46GAixEiYiI8mLz5s0YP3585mOOiBKRLixEySrs3bsX7du3BwBUqFABERERmhMREVmmtLQ0TuhGWUpJSUFgYCBq166duaxChQoaExFRUcZClKxC69at0bdvXwDAqlWr0KNHD82JiIgsU2xsLBYvXqw7BlkgBwcH+Pj4ICoqCr1790br1q3h4OCgOxYRFVF2ugMQmSo8PBwAUK9ePZQtW1ZzGiIiy5SUlITPP/8cnTp1gouLi+44ZGH8/PwQFxeHkSNH4qWXXtIdh4iKMI6IktUIDw9H5cqVWYQSEWXj7t27iImJwZdffqk7Clmg+fPno0GDBnB3d0fx4sV1xyEqMOfPn8eJEyd0x6BscESUrEZ4eDheffVV3THIQt25cwc//vgjSpUqBVdXV9SqVQuurq4oVqyY7mhEZpWUlAQAWLBgAQYNGoQ33nhDcyKyFEePHsWhQ4ewcOFC3VGKhISEBKSlpcHOjj+3dQgICAAAvPzyy5qT0JPkOCKqlKqulNqjlDqplDqulPokizZKKfWDUuqsUuqYUur13IQ4dOgQBg4ciP/+97/YsWMHoqKicvN0KgJEBMeOHWMhqklqairmz5+PzZs349q1a7rjZKlEiRLo3Lkz/vvf/6J9+/Z4+eWXUbx4cVSuXBmDBw9Gamqq7oj0ABHBlClT4Ovri7Nnz0JEdEcqNO4XoiKCYcOGwWAwaE5ElmL+/PkoVaoUevbsqTtKkfDPP/+gUqVKGDhwINavX5/ZN8k8AgICsHXrVt0xKBumHJqbBmCUiLwEoDEAb6XUo7sW2gKok3EbAmBubkI0btwYrVq1wueffw5PT09UrFgRlSpV4mQLBcRgMCAoKAhxcXG6o5jsypUriIuLYyGqib29PZo1a4ZBgwahSpUqqFy5Mtq3b499+/bpjvaQF154AYcOHcKbb76Zuez69et44YUXWIhqJCJITExEZGQkTp06heDgYOzatQtXrlzBoEGDUKdOHVStWhU9evTAnDlzEB4ezllfn8KDP3bDwsKwYMECjWkKt4MHD1pNcZGYmIjly5fj/fffR8mSJXXHMVlQUJDuCHlWp04dvPfee9i+fTs6d+6McuXKoVOnTvjll18QHR2tNVtycjIOHTqkNUNBunHjBoKCgrB//34kJCTojpNraWlpuiOYxGg04sqVK3l+fo7HCojINQDXMu4nKKVOAqgK4MGDrjsCWCrpu7QPKaWeUUpVzniuSXr37o3ixYujR48eSElJQXx8PFq0aJG7f41G8fHxKFWqlO4YT3Tjxg1s27YNW7duxe+//w6j0Yi///4bpUuX1h3NJMePHwcAqypEjUYjjEYjDAYDHB0ddcd5aq+88gr27t0LDw8PXLt2DQcPHsTrr+fq4AezKF++PHbt2oUBAwZg9erVqFq1KsaMGYMZM2bgjz/+QJ06dXRHLLS++uorhISEID4+/qFbQkJCjoXltWvXsHbtWiQmJsLFxQW1a9cu8MOqt2zZgjt37qBNmzZwdnYu0Ncyp6SkJHTv3h1+fn6YNWsWPD09dUcyidFohI2NdU1dMWHCBAQFBaFNmzbo1KkTvLy8LHaCqC1btiAhIQFDhw7VHcVkBw4cQJ8+fTBp0iQMGDBAd5xcK1WqFObOnYs5c+YgJCQE69evx4YNG/DBBx/AxsYGXbp0wZo1a8yWx2g0IigoCMuXL4efnx9Gjx6N11577bHP7CfdSpcu/dA13S3Z7t278dxzz0EphcDAQLRr1053JJMsW7YMHh4eCA8Px4kTJ/Dvf/8bSindsZ4oKCgI/fv3x549e1CjRo3cb0BETL4BcAVwCUCpR5ZvAtDsgce7ADTI4vlDAIQACKlRo4ZkZdu2bVKsWDEZOXKkGI3GLNtYGqPRKHXr1pXhw4frjvJEFy9elFmzZknz5s3FxsZGypQpI2fOnNEdy2RGo1GuXLkiKSkpuqPk6N69exIRESEuLi4ycuRIGTVq1BPbAgiRXPTBgr6Z0kdPnz4tVatWlenTp+fp/TEXg8Eg48aNk59//llCQ0Nl1KhRYjAYdMcySWhoqLRr106aNGkiixcv1h3HZMOHD5cWLVpIx44dpW/fvuLt7S1ffPGFTJs2TebMmSPLli2TDRs2yJ49eyQ0NFRmzZolAOSNN96Q77//Xm7cuGHWvF26dBEA4uDgIG3btpX58+dLZGTkQ22ssY/GxMSI0WiUM2fOWM33qIhInTp1ZPz48bpj5Mobb7whAKRKlSoyaNAgWb9+vcV+zhiNRjl69KjuGCYzGo3y+uuvCwABIHPmzMmynTX20dOnT8vMmTPlm2++yfP7kxsnT56UcePGiaura+b7CUBsbGweepzVTSklpUuXlurVq4unp6dZ8uYHo9EoiYmJcufOHUlLS9MdxyQxMTHi5OQkdnZ20qRJEwEgw4YNk9TUVN3RspSWlia1atUSAOLq6ioXLlzIsl12fTQ3HaskgFAAXbJYtzmLQtQ9u+25u7s/8R+2b98+uXz5ct7fGTM7ePCgAJBVq1bpjmKS6OhoWbx4sSQkJOiOUigFBQVJ9erVMz/Evb29n9jW0r5AH7xl10fPnj0rsbGxuXtjNImLi9MdIdeSk5PF2dlZAMiBAwd0xykwmzdvllOnTml7faPRKGFhYeLj4yP16tXL7LONGzeWadOmycmTJ622j1qbuLg4ASCzZs3SHSVXli1bJkePHrWqgt9ahIeHS5UqVR4qimbMmPFYO/bR7CUnJ8v69etl1KhR8sYbb4itrW3m+9mmTRv54YcfZPHixbJu3TrZuXOnBAcHy6lTp+Tq1auSkJDAv20zi4yMlHHjxj32/8kSf8scPnxY6tatK0opASDVq1eXs2fPPtYuuz6q0tdnTyllnzHquV1E/pfF+vkA9orIqozHpwE0l2wOzW3QoIGEhITk+NrWoF+/fti6dSuuXLlSKA7BpKc3ePBg/PzzzwCA0aNHY+bMmVm2U0qFikgDc2YzVWHqo9aoS5cuOHjwICIjI2Fra6s7TpFw4cIF/Pbbb1i/fj32799//3Bi9lEzCAoKQtOmTbFr1y54eHjojkMW4t69ezh58iSOHTuWeevYsSO8vb0zD1fk92juJCQk4I8//kBgYCCuXLmCRYsW8TvGgogIRo4ciVmzZj203M3NDZs2bUL16tU1JXuyhIQEHD16FGFhYYiKisIXX3yBEiVKZK7Pro/meI6oSu/pPwM4mVURmuE3AMOVUr8CaAQgLrsitDC5efMm/Pz8MGLECBahlGn69OkICAjAzZs3efkQypO2bduiXLly/IFgRrVq1cInn3yCTz75BLGxsdi8eTP69++vO1aREB4eDsC65gGggufk5IT69eujfv36mctMGUChJ3N2doanp6fVnDte1Cil8N1332Hy5Mm4ePEiLly4kPnf6dOnY/LkyShbtqzumA9xdnZGs2bN0KxZs1w/15QLGzUF0BdAuFLqSMayLwHUAAARmQdgC4B2AM4CSAIwMNdJrNTSpUuRnJyMIUOG6I5CFqRcuXKYNm0ahg4dyguGU560bdsW1apV0x2jyCpbtiz69evHQtRMwsPDUbFiRZQvX153FLJwljxxC1F+cXZ2xquvvlrod86ZMmvuAQDZ9vqM43+98yuUNTAajVBKYf78+WjZsiVq166tOxJZmPuH57IQpbyoVq0aqlSpojsGkVmEh4cX+h9cRET0MOuaJ92CTJo0Cb/++itOnTplVVOhk/nY2Njgp59+KlSXhSDzsrZLWRDl1qVLl5Camorw8HC4ubnpjkNERGbEXzl59Ndff6F3796ws7PD1atXce7cOd2RyAK5u7uja9euumMQEVmkY8eOoWHDhrh16xZu376NRYsW6Y5ERERmwkI0j27cuAEASEtLw65du1CrVi3NichScUSUiChrNWrUwJEj6dNP+Pr68rxoIqIihIVoHl2/fh0A8MILL2DZsmU8hI6IiCiXatSokXm/UaNGaN26tcY0RERkTqye8sBoNCIqKgrOzs5Yv349SpcurTsSERGR1SldunTmUSOTJk3ijKhEREWIKZdvoUfcvHkTaWlpWL58OV588UXdcYiIiKySUgo1atTIvLYhEREVHSxE8+DGjRuYMGEC3n33Xd1RiIiIrFrNmjXxySefcDSUiKiI4aG5eVClShX4+PjojkFERGT1evbsiVatWumOQUREZsYR0TwoU6aM7ghERESFwvvvv8/RUCKiIogjokRERKQNi1AioqKJhSgRERERERGZFQtRIiIiIiIiMisWokRERERERGRWLESJiIiIiIjIrFiIEhERERERkVmxECUiIiIiIiKzyrEQVUr5KqWilFJ/P2F9c6VUnFLqSMbNJ/9jEhERERERUWFhZ0KbxQBmA1iaTZtAEfHKl0RERERERERUqOU4Iioi+wHcNEMWIiIiIiIiKgLy6xzRJkqpo0qprUqpV57USCk1RCkVopQKiY6OzqeXJqL8wj5KZNnYR4ksG/sokenyoxANA1BTROoB+BHA+ic1FJEFItJARBqUL18+H16aiPIT+yiRZWMfJbJs7KNEpnvqQlRE4kUkMeP+FgD2SqlyT52MiIiIiIiICqWnLkSVUpWUUirjfsOMbcY+7XaJiIiIiIiocMpx1lyl1CoAzQGUU0pdATABgD0AiMg8AN0ADFNKpQG4C6CniEiBJSYiIiIiIiKrlmMhKiK9clg/G+mXdyEiIiIiIiLKUX7NmktERERERGbAgw+pMGAhSkRERERkBRITEzF27FjEx8frjkL01FiIEhERERFZuL1798LNzQ23bt1C6dKldcchemosRImIiIiILNSdO3fw8ccfo0WLFrhw4QK8vb11RyLKFzlOVkREREREROa3b98+fPDBBzh//jwA4O2334abm5vmVET5gyOiREREREQWRkRgNBofOgyXo6FUmLAQJSIiIiKyMEopJCQk4MiRI3B3d0flypXRuXNn3bGI8o1VFqIiAl9fX6SmpuqOQkRERESU78LDw9GnTx+8/fbbOHDgAGbOnAl7e3vdsYjyjdUVoomJiejduzd27tzJzkhEREREhU50dDTeffddlC9fHv7+/nBycsL777+vOxZRvrKqQvSff/5B48aN8euvv6J79+664xAREVmlI0eO6I5ARE+QkpKCbt26ISYmBr/99hvKlSunO1KW/vrrL90RyMpZTSEaEBCABg0a4Pjx4yhVqhTatGmjOxIREZFVSUpKgre3N9auXas7ChFlQUQwfPhwBAYGYuXKlahbt67uSFny9/fHkiVLdMcgK2cVheiqVavQt29fJCQkAAA6duwIJycnzamIiIisx+HDh+Hu7o6ffvoJAwcO1B2HNElLS9MdgbIxe/ZsLFy4ENOmTUOHDh10x8nSkSNH0L9/f1SpUkV3FLJyVlGI9urVCz169ICNTXpcHpZLRGRd0tLScPjwYcydOxcDBw5EcHCw7khFhtFoxMyZM9GoUSOcOnUK77zzDp599lndsUiT77//HgaDQXcMysKOHTvw6aef4v3338eYMWN0x8lSVFQUOnbsiKSkJFStWlV3HLJydroDmMLf3x++vr6YMWMGQkND4enpqTsSET2BiCA2NtZiz2kh8zl8+DD8/Pzwxx9/4K+//kJSUhKUUli8eDEaNmyoO16RcO/ePfTu3RsBAQGZywYPHqwxEekUEREBHx8feHt7w9bWVnecPElKSsKdO3dQvnx53VHy1T///IPu3bvjjTfewMKFC6GU0h3pMSkpKejatSsuXboEABwRpadm8SOiV65cwZAhQ9C8eXOMGjUKixcvhoODg+5YlM+Cg4MRFBSEiIgIpKSk6I5DuSQiCAsLwxdffIHXX38dN27c0B2JLMALL7yAQ4cOYd++fUhKSgIAzJ8/H/369dOcrOhwcnLC3LlzM3cMubi4oFOnTppTkQ73zz1MSkpCcnKy7ji5FhYWho8++gitWrVCyZIldcfJV7dv30aHDh1QokQJBAQEWOTpZyICb29vHDhwIHMZR0TpaeU4IqqU8gXgBSBKRB47Y1ql77L5HkA7AEkABohIWH6EMxqN6N+/P0QES5cuha2tLYoVK5YfmyYL4+rqim7duiEwMBAAUL58eVStWhVVq1aFj48PR08skIggNDQUa9asgb+/P86fPw8A8PPzwyuvvKI5Hel24sQJjBkzBnv37s1c9uOPP+LDDz/UF6oISktLQ69evZCcnIzvvvsOFy9etMgfuVTwAgICsGnTJgCwmh2+t2/fxsqVK7Fo0SIcPnwY9vb2CA4OLlS/BdPS0tCzZ09cvnwZgYGBqFy5su5IWRIRTJ8+HdeuXUNwcDDi4uJYiNJTM+XQ3MUAZgNY+oT1bQHUybg1AjA3479P7X//+x92794NPz8/VK9ePT82SRaqQoUK2LlzJ4YPH46FCxciOjoa0dHRsLe3R40aNXTHoyxERkZi6tSpDx3yN3r0aLz33nsaU2Xv/PnzKFeuHEqVKqU7SqEVFRWFCRMmYOHChShbtizmzZuHwMBA1K9fH8OHD9cd7zGrVq2Ch4cHKlasqDtKgfjqq6+wZ88erF27Fp07d8b169d1RyIN4uPj8fHHH2c+toYR0RMnTqB9+/a4ePFi5rKpU6fitdde0xeqAIwePRrbt2+Hn58f3N3ddcd5IhsbGyQmJmLr1q2YMmUKXnzxRX6X0tMTkRxvAFwB/P2EdfMB9Hrg8WkAlXPapru7u2QnLCxM7O3tZcCAAdm2o8LFaDTKjz/+KLa2tpk3e3t76dmzp4SFhemOl+8AhIgJfVDHLac+eu7cOWnUqJEAEADi4eEhqampT/2eFKTt27dLlSpVZO3atWI0GnXHyVFCQoLuCCZLS0uTqVOnirOzszg5OcmXX34pcXFxIiJy8uRJzemydvHiRRkyZIjY2dlJly5dZOvWrZKWlvZQG2vuo2vXrhUA8vnnn+ftDTKz27dvS2Jiou4YhdI333wjdevWzfy8PnfunO5IOUpOTpZ+/fplZm7RooUYDIbH2llzH120aJEAEB8fn9y/QRqMHTtWHBwcJCoqSncUk507d+6xz/X7wsPDLeq3wN27dzO/Ny1dWlqaye9ddn00P84RrQrg8gOPr2Qse4xSaohSKkQpFRIdHZ3tRp2cnODl5YUffvghHyIWvLNnz6JPnz6Ii4vTHcUkQUFB8PHx0R3jMUopDB8+HNu3b8fzzz+Pixcv4vPPP8eePXtw4sQJ3fFMtnjxYnz99dcwGo26o+RKbvqoUgo3b97EN998gxo1auDXX3+FnZ3lzn+WmpqKrl27IjIyEl27dkXHjh0zJ1ywREajEe7u7vDw8MDq1ast/lA6Gxsb7Ny5Ex07dsTp06cxZcqUzL3lL774ouZ0WRs/fjwWLFiAtLQ0rFu3Dm3btsWzzz6LSZMm4fLlyzlvQIPc9NEyZcqgR48e+Prrr82U7ul4e3tjypQpumOYbNSoUdi9e7fuGCb57LPPsHz5cowZMwbffvstUlNTdUfKka2tLWJiYjBo0CA888wzWLJkSebVEyxZbvpo2bJl0b9/f0yYMMFM6Z5Oo0aNMHHiRKuZKCohIQHPPfccnJ2d0aBBAwwcOBD/+9//sGPHDly/fh2+vr5o0qQJtm3bdn8wTatJkyZZzWRyYWFh6N27d+YpdXn2pAr1wRuyHxHdDKDZA493AXDPaZs57SWyNidOnBBbW1sZO3as7igmmTVrlgCQe/fu6Y7yRBcvXsy8n5ycLCkpKRrTmO7GjRtSpkwZ8fLyynFvEax4T65I+h6xpKQk+euvv3J+YzS7d++e1KtXL3PvOgApUaKEfPvttxY5kpuSkiI//vijvPLKKwJAypcvL2PGjJGzZ8+KiMjOnTslPDxcc8qHJScn646QK8eOHZM333wz8++hQoUK4uXlJZMnT5YdO3aIwWCw+j5qTW7fvi137tzRHSNHwcHBEhoaKgBk3rx5uuPkiSWNAmXHYDBISkqKrFu37olt2EfpSZKSksTPz0/Gjx8vnTp1kueee06UUpmf+Y6Ojpn3GzduLNu2bdPaN+7duye3b9/W9vq58eOPPwoAqVu3rly/fj3bttn1UZW+PntKKVcAmyTryYrmA9grIqsyHp8G0FxErmW3zQYNGkhISEiOr21NPvroI/j6+uL06dOoWbOm7jjZ+uWXX/DBBx8gOjqal9nIZ3379sW6detw4sSJHP8OlFKhItLATNFypbD10TNnzsDPzw8ODg5wcHCAo6Nj5v3XX38dL7/8su6IWRIRBAUFYd68eVizZg2Sk5PRsmVL1KlTB0uXLsXixYvRrVs33TGt0s2bNzFz5ky4u7ujYcOGqFGjxmOXTGAfUFXjsQAAGQtJREFUpUe9+OKLcHZ2xoULF3Dp0iUUL15cd6QijX2UcuPOnTs4fvw4wsPD4ePjg8jIyIfWN2nSBP/5z3/wzjvvaEpoHfr27YsdO3Zg9+7dOf5+yq6P5sdxdL8BGK6U+hXpkxTF5VSEFlYTJ07E8uXLMW7cOCxfvlx3nGzdP2QuPj6ehWg+2rVrF5YvX46ZM2da/M6IoqZOnToYN26c7hi5ppRC06ZN0bRpU8yaNQtLly7F/PnzsXPnTgDAe++9hy+++AL/+c9/rPa6gLq4uLhg+vTpumOQFTEYDDh37hzS0tLg7OyMrl27wt/fHyVKlNAdjYhMUKJECTRs2BBGoxF169ZF+/bt/1979x5VVZ32Afz7cDhKoYhmSoJg3l5lqVQqps44jdiEmlqaUypY+ppT2ShTJpqay9Zq8F1qqV2WZReZN2VyjZc0MaVSYxIUvGVortFKw2te8QaCPO8fKC8ZwhHPOb+9z/l+1jrLsw+HwxfkOZtn//b+/RAZGfmrW0REhOmYlpefn48NGzbc8qU3rizfkgbgAQANRSQfwDQATgBQ1fkA0lG2dMs+lC3fMuKWEtlYo0aNMGnSJLz88ssYN24cIiIiLDsNd8VGlNyjsLAQzz77LNq3b4+kpCTTccgH3XHHHfjb3/6GZs2aYfDgwbhy5QoAICUlBdu3b8fixYtRv359wymJfNehQ4dQUlICoOwaxjfeeINNKJEN3X///Vi7dq3pGLZUUlKCBQsWoGXLlrf8WtU2oqo6pJqPK4Axt5zERyQlJWH+/Pl48cUXoar45JNP0KRJE9OxfoONqPvNmDED+/btwzfffAOn02k6DvkoVUWjRo2QlpaGo0eP4tixY+X/Pvfcc3jjjTcQFhZmOiaRT/rxxx8BAE6nE8uXL7fsRFxERJ4SGBjoliYUcM+puVRBZmYmOnXqhGXLlgEAsrOzMXDgQMOpfouNqHvt3bsXKSkpGD16NLp27Wo6Dvmwa6fqEpH3XWtEFyxYgAceeMBsGCIim7P+PNg2ExERgU2bNpVvZ2dnG0xzY2xE3UdV8dxzzyE0NBQpKSmm4xARkYf8+OOPmDp1Kp588knTUYiIbI8jom4WHR2NjRs3omfPnjh06BCysrJMR6oUG1H3+fjjj/HVV19h0aJFvD6PiMiHde3aFQ899JDpGEREPoEjoh7QunVrfP3114iKikJubq4lF6IPDg6GiLARraFr6x+dOnUKL774Inr16oUhQ6q8nJqIiGwuPj7+N0v8EBFRzbAR9ZDmzZtj48aNCA8Px7fffms6zm8EBASgbt26bERraP/+/Vi+fDmSk5NRUFCAd955h3+cEBERERG5iKfmelBUVBQ2btyI3bt3m45SqZCQEDaiNZSRkYGpU6fi5MmTmD59Olq1amU6EhERERGRbXBE1MPCw8Px4IMPmo7xK8ePH0dubi5CQkLwyy+/YN68eeXrEZJrMjIycPLkSQBAWloaNmzYYDYQEREREZGNcETUDzVo0AAxMTE4fvw4du/ejfz8fIwdO9Z0LNsoKSnBV199BaBsLaWXX36Z0/gTEREREd0Ejoj6ocDAQAwZMgSlpaUAgB49ehhOZC85OTk4e/Ys6tati/T0dCQmJpqORERERERkK2xE/VRCQkL5fTaiNycjIwNNmjRBZmam5U67JiIiIiKyAzaifuree+9FdHQ0AgIC0K1bN9NxbOXIkSPIyspCTEyM6ShERERERLbERtRPiQgSEhJw7733om7duqbj2MrMmTMRGRlpOgYRERERkW1xsiI/NnToUJw4ccJ0DNupU6eO6QhERERERLbGEVE/FhUVhaSkJNMxiIiIiIjIz7AR9XNNmzY1HYGIiIiIiPwMG1EiIiIiIiLyKpcaURGJF5G9IrJPRCZW8vGnROQXEdlx9TbK/VGJiIiIiIjIF1Q7WZGIOAC8DeBBAPkAckRkparuvu6pn6jq8x7ISERERERERD7ElRHRWAD7VPUHVb0M4J8ABng2FhEREREREfkqVxrRcAA/V9jOv/rY9QaJyLci8i8RqXQGHBEZLSK5IpL7yy+/1CAuEXkSa5TI2lijRNbGGiVynSuNqFTymF63vQpAM1XtAOALAKmVvZCqvqeqnVS105133nlzSYnI41ijRNbGGiWyNtYoketcaUTzAVQc4YwAcLjiE1T1pKoWXd1cAKCje+IRERERERGRr3GlEc0B0EpE7haRWgCeALCy4hNE5K4Km/0B7HFfRCIiIiIiIvIl1c6aq6olIvI8gLUAHAA+VNU8EXkVQK6qrgQwVkT6AygBcArAUx7MTERERERERDZWbSMKAKqaDiD9usdeqXB/EoBJ7o1GREREREREvsiVU3OJiIiIiIiI3IaNKBEREREREXkVG1EiIiIiIiLyKjaiRERERERE5FVsRMnjCgsLkZ6ejpKSEtNRiIiIiIjIAtiIkscFBQVh27ZtaNGiBWbMmIETJ06YjkRERERERAaxESWvSE5ORsOGDTFp0iRERETgqaeeQm5urulYVEOXLl3CpUuXTMcgIiIiIptiI0pe4XQ6sXDhQjidThQVFSE1NRUDBw7El19+aToa1cDWrVuxbNky0zGIiIiIyKbYiJLXtG/fHtOmTSvfbty4Mdq0aWMwEdVUdnY2PvjgA9MxiIiIyEUZGRn44YcfTMcgKsdGlLxqwoQJuO+++zB58mTs378fnTt3xpYtW0zHopuUnZ2N9evXc4dGRERkE1FRUYiNjUWLFi3wzDPPYOnSpTh9+rTpWOTHLNWIFhUV4eeffzYdgzzI6XQiNTUVEyZMwJYtWxAaGooePXpg8eLFpqPRTcjOzgYAfPTRR4aTkCeVlJSgsLDQdAxy0bp163D06FHTMYjIolq3bo1Vq1bh8OHDePfdd/HYY4+hYcOG6N+/P86ePWs6HvkhSzWigYGBSEhI4DIfPq5du3YICQlBy5YtkZWVhbi4OAwbNgyTJ09GaWmp6XhUjfz8fBw6dAgAsHDhQly5csVwIvIUh8OBiRMnYvjw4fjss89w+fJl05GoCrVq1UJ8fDzOnDljOgoRWVTXrl2RlpYGEQEAlJaWIi4uDiEhIYaT+TZV5SSPlbBUI+pwOLBt2zZMnz7ddBTyknr16mHlypUYP348/v73v2PgwIE4d+6c6VhUhaysrPL7+fn5yMjIMJiGPElEkJKSgl27dqFfv35o3LgxRowYgc8//xzFxcWm49F12rZti507d6J///78g4ewfv368oOGZE2mDuQ+8sgjePPNNwGUDQ4kJSUhLi4Oe/fuNZLHH4gIZs2ahdGjR2PLli1QVdORLMFSjShQ1pi89tpr+OKLL0xHIS9xOByYOXMmFi5ciDVr1qB79+746aefTMeiGzh79ixSUlIAoHwmZPJdt912G5YuXYrQ0FCcOXMGCxcuxMCBAzF37lzuSC2mUaNGaNCgATIzM/H444/zYIGfi4iIQMeOHfH111+bjkI38P777xv72mPGjMGECROwevVqLF68GHl5eejQoQOmTZvGSzI8ZMKECcjJyUGXLl0QExODefPm4dSpU6ZjGWW5RjQkJASqioSEBBw7dsx0HPKiJ598EuvXr8exY8fQuXNnZGZmmo5ElRg1ahQGDRoEoKxe4+LiDCciT2vevDk+/vjj8u3AwECEhYWVn9pF1iAiaNu2LQBg1apVGDVqFC938GOtWrVCWFgY4uLiMG/ePB44spgjR45g+vTpRi9HS0lJQXh4OIYMGYLvv/8eI0aMwKuvvooOHTpwQMgDateujUWLFiEoKAi7du3CuHHj0KRJE6SkpPhtfVqyEQWAY8eOISEhgTtRP9OtWzfk5OQgIiICcXFxXCLEoho2bAgAOHHihOEk5C19+/bF1KlT0aVLF3Tr1g2JiYkYNmwYr0e0mGuNKAAUFxeXTyxG/mnw4MEoKSnBuHHjkJCQgIsXL5qORFdt3rwZR44cweeff24sQ0BAABwOBwCgfv36mD9/PjZt2oSgoCA8+OCDHBTygOjoaMycObN8u2HDhkhMTPTbA7suNaIiEi8ie0Vkn4hMrOTjtUXkk6sf3ywizWoaKCQkBLVr10atWrUwdOhQHDx4sKYvRTYVGRmJf//73+jfvz9GjRqFpKQkTmBlMfXq1UNAQABOnjxpOgp50bRp0zBlyhSkp6dj7ty5WLp0Ke655x6evWAh0dHRGDZsGEQEMTEx6Natm+lIN1RSUoLt27dzEiwPGjx4cPn9xYsXIy4ujst1WMS1g0RWO+DetWtXbN26FTNnzsTy5cvRpk0bvPfee7YbGFJVy06mOGbMGPTp0wcdOnTAxYsX0bFjR2zYsMF0LCOqbURFxAHgbQC9AUQDGCIi0dc97b8BnFbVlgDeAPA/NQ3Ur18/LFmyBJcvX0ZYWBiaNWtW05ciGwsODsaSJUvwyiuvYO7cuejbt68tRl5UFTk5OXj22WeRn59vOo7HBAQE4I477uCIqJ9xOBx4+OGHERAQgLFjxyInJwchISF44IEHMGXKFFtck/jdd9+ZjuBRgwYNQmpqKgYNGoTZs2dbegQsMDAQ3333HUJDQ9G1a1ckJSUhLy/PdCyf0rp1a7Rv3x4AcPvttyMtLQ3169c3nOr/FRcXY+XKlXjllVdMR/G6a43oZ599ZrlRR6fTifHjx2P37t34/e9/j7/85S/43e9+h127dpmO5pI9e/agd+/elj3IJSL48MMP8cILL2Dr1q0IDw9Hr169MHv2bL87RdeVEdFYAPtU9QdVvQzgnwAGXPecAQBSr97/F4A4qeEY81//+lf07dsX4eHh2Lp1a01ewoiTJ0+ioKDAdAyXHThwABcuXDAdo0oBAQGYPn06lixZgszMTEvvqE6cOIE5c+agQ4cOiI2NRWpqKrZv3246lkdFRUXZ6g1TVbF69Wr89NNPWLNmjek4VXr77bexe/duFBUVmY5Spfbt22PLli0YO3YsXnvttfJZGK2otLQUM2bMwMiRI7Fjxw6MHz/eVvsYV0VGRsLhcGDy5MkoKCj41SzXVpSYmIh33nkHmzdvxvvvv2+bg8+LFi3Cnj17bHG2zuDBg5GWlobbbrsN//jHP0zHAQDs3bsXycnJiIyMxIABA7Bz507MmjULkyZNwtNPP42BAweiR48ePlmjQNnZANcGW8aMGYPdu3ebjlSpqKgofPrpp1i2bBkOHjyIbt26WXpQ4MKFC5g4cSJiYmKQlZWFMWPGmI50Q40bN8bw4cNx991345tvvsHw4cMxfvx4S69EUFpaijlz5gAoWzXhpZdeuvWVLlS1yhuAxwC8X2E7EcBb1z3nOwARFbb3A2hYyWuNBpALIDcyMlKrcv78+So/bhXFxcX67rvv6p///GeNiYnRgwcP6unTp03Hqlbfvn21U6dOpmO4bMeOHVpQUGA6RqUOHDigo0eP1jZt2igABaAioqtXr672cwHkajU16M3bzdSo3axatUoBaJs2bbR///6m41SpdevW6nQ69fHHH9e8vDzTcVyyfv16LSwsNB2jUseOHdM//elPCkADAgIUgN5+++360UcfVfu5dq7RkydPuvLjsYQFCxbo/PnzTcdwWffu3dXhcOgzzzyj+/fvNx2nStf2nQcOHNDS0lJjOa5cuaKLFi3S7t27l+8rr7+FhoZqy5YttUuXLtq3b1/Nzc2t9nXtWKMm/x9qqqCgQNetW2c6RqVKS0t12bJlGhkZ+avfp5iYGNv0E6WlpbphwwbTMW7o8uXLOmzYMG3cuLGOHj1anU6nBgUFaUZGRrWfW1WNulJQgytpRN+87jl5lTSid1T1uh07dnTTj8a89u3bKwANCgrSpk2banFxselIVTp37pzWrl1bp0yZYjqKzzl16pSuWbNGp02b5tIbitV2oBVvvlSj586d09dff13r1q2rANThcOihQ4dMx6rU+fPnVUTK31OsvGOygy+//FLDwsJ+9cdJYmKinjt3zqXPZ416j9X3ndeUlJRocHCwAtDY2Fg9c+aM6Ui2UVhYqNnZ2Tpnzhx94okntFmzZgpAe/ToUeP/f9aofzt//ry+8MILGhkZqXfeeacGBweX70P79Oljm/cVK7tw4YL27t27fB8aHBysL730kh49etSlz6+qRgNdGDTNB9C0wnYEgMM3eE6+iAQCqAfAbxbGiY2Nxa5du1BYWIjhw4cjMNCVH6s5a9euRVFRER555BHTUXxO/fr1ER8fj/j4eNNRqII6deqgefPmqFWrFoCyRcRTU1MxadIkw8l+Ky8vD6oKp9OJ5cuX4w9/+IPpSLZ05coVzJo1C3PnzkXdunURHh6OkJAQ1KtXD06nExcvXkSdOnVMx6QKrL7vvOY///kPLly4gI4dO2Lt2rWoV6+e6Ui2Ubt2bXTp0gVdunTBuHHjAABHjx7F5s2bUVBQgAYNGhhOSHYTHByM2bNnY/bs2eWPqSqKiopw6dKlawNkVEOnTp1Cv379sGnTpvLHQkND8eijj6Jx48a3/PquvOvnAGglIncDOATgCQBDr3vOSgBPAshC2am8X6kf/c/HxsaWz3o2cuRIw2mqt2LFCkREROC+++4zHYXIawYMGIDY2FiMGDECa9euxQcffIDk5GQEBFhrFatvv/0WDocDn3zyCQ9o3AKHw4Hk5GQkJyebjkI+Ztu2bbjnnnuwbt06hIaGmo5je2FhYRgw4PqpR4hqTkQQFBSEoKAg01Fs7fDhw4iPj0dBQQF69+6Ntm3bIjo6Gm3btkWbNm3c8jWqbURVtUREngewFoADwIeqmicir6JsqHUlgA8A/K+I7EPZSOgTbklnE507dwYA9OrVC82bNzecpmrFxcVYvXo1hg4d6rdrFpH/uuuuu5Ceno633noLEyZMwMaNG/HHP/7RdKxf2bVrF1JTU/Hoo4+ajkJElSgpKUFGRgZH74jIp4kIsrKyEBwc7LGv4dJ5MKqaDiD9usdeqXC/EGXXkvqldu3aISgoCE8//bTpKNXKzMzE6dOnefSR/Na1pUd69uyJ9PR0yzWiI0eORExMjOkYRHQD/rz4PBH5j7vuusvjX8MeF2RYnNPpxEMPPWSL5u7TTz9FvXr1eN0Z+b127dohOvr6JZHNYxNKZG1sQomI3MNaF0fZ2OzZs1G7dm3TMaqkqlixYgX69OlTPmkLkT+z2vWhRERERP6Cf4W5SYsWLUxHqNbOnTtx8OBBzpZLRERERERGsRH1IytWrIDT6eRMnEREREREZJSYWmVFRH4BcKCapzUEcMILcdyFeT3LF/NGqeqd3ghzs1ijlsC8nsUatR7m9SxfzMsa9S7m9SxfzHvDGjXWiLpCRHJVtZPpHK5iXs9iXuux2/fIvJ7FvNZjt++ReT2Lea3Hbt8j83qWv+XlqblERERERETkVWxEiYiIiIiIyKus3oi+ZzrATWJez2Je67Hb98i8nsW81mO375F5PYt5rcdu3yPzepZf5bX0NaJERERERETke6w+IkpEREREREQ+ho0oEREREREReZVlG1ERiReRvSKyT0Qmms5TFRH5UESOi8h3prO4QkSaish6EdkjInkiMs50pqqISJCIbBGRnVfzTjedyRUi4hCR7SLymeks7man+gRYo57GGrUe1qhnsUa9gzVqHaxRz7JjjbqjPi3ZiIqIA8DbAHoDiAYwRESizaaq0kIA8aZD3IQSAC+qalsA9wMYY/GfbxGAnqoaA+AeAPEicr/hTK4YB2CP6RDuZsP6BFijnsYatRDWqFewRr2DNWodC8Ea9SQ71ugt16clG1EAsQD2qeoPqnoZwD8BDDCc6YZU9WsAp0zncJWqHlHVbVfvn0PZL1G42VQ3pmXOX910Xr1ZepYtEYkA0BfA+6azeICt6hNgjXoaa9RyWKMexhr1PNaotbBGPctuNequ+rRqIxoO4OcK2/mw8C+PnYlIMwD3AthsNknVrg7/7wBwHECGqlo6L4A5ACYAKDUdxANYn17EGvUY1ii5BWvUY1ij5BasUY9wS31atRGVSh6z7FEBuxKROgCWAkhS1QLTeaqiqldU9R4AEQBiRaSd6Uw3IiIPAziuqltNZ/EQ1qeXsEY9gzVK7sIa9QzWKLkLa9T93FmfVm1E8wE0rbAdAeCwoSw+SUScKCvMRaq6zHQeV6nqGQAbYO3rFLoD6C8iP6HsdJueIvKx2Uhuxfr0AtaoR7FG6ZaxRj2KNUq3jDXqMW6rT6s2ojkAWonI3SJSC8ATAFYazuQzREQAfABgj6q+bjpPdUTkThEJvXr/NgC9AHxvNtWNqeokVY1Q1WYo+939SlUTDMdyJ9anh7FGPYs1SreKNepZrFG6VaxRz3FnfVqyEVXVEgDPA1iLsouLl6hqntlUNyYiaQCyAPyXiOSLyH+bzlSN7gASUXYEY8fVWx/ToapwF4D1IvItyt68M1TV56Zytwu71SfAGvUC1qiFsEa9gjVKNcYa9QrWqA2IKk9JJyIiIiIiIu+x5IgoERERERER+S42okRERERERORVbESJiIiIiIjIq9iIEhERERERkVexESUiIiIiIiKvYiNKREREREREXsVGlIiIiIiIiLzq/wAj3L7bYZ97HgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1152x288 with 4 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "\n",
    "feat_type_in = nn.FieldType(r2_act, [r2_act.trivial_repr])\n",
    "feat_type_hid = nn.FieldType(r2_act, 8*[r2_act.regular_repr])\n",
    "\n",
    "model = nn.SequentialModule(\n",
    "    nn.R2Conv(feat_type_in, feat_type_hid, kernel_size=3),\n",
    "    nn.InnerBatchNorm(feat_type_hid),\n",
    "    nn.ReLU(feat_type_hid),\n",
    "    nn.R2Conv(feat_type_hid, feat_type_hid, kernel_size=3),\n",
    "    nn.InnerBatchNorm(feat_type_hid),\n",
    "    nn.ReLU(feat_type_hid),\n",
    "    nn.R2Conv(feat_type_hid, feat_type_out, kernel_size=3),\n",
    ").eval()\n",
    "\n",
    "S = 11\n",
    "x = torch.randn(1, 1, S, S)\n",
    "x = nn.GeometricTensor(x, feat_type_in)\n",
    "\n",
    "fig, axs = plt.subplots(1, r2_act.fibergroup.order(), sharex=True, sharey=True, figsize=(16, 4))\n",
    "\n",
    "X, Y = np.meshgrid(range(S-6), range(S-7, -1, -1))\n",
    "\n",
    "# for each group element\n",
    "for i, g in enumerate(r2_act.testing_elements):\n",
    "    # transform the input\n",
    "    x_transformed = x.transform(g)\n",
    "    \n",
    "    y = model(x_transformed)\n",
    "    y = y.tensor.detach().numpy().squeeze()\n",
    "    \n",
    "    # plot the output vector field\n",
    "    axs[i].quiver(X, Y, y[0, ...], y[1, ...], units='xy')\n",
    "    axs[i].set_title(g*90)\n",
    "    \n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We want to mention that vector fields can not be acted on by ReLU nonlinearities since this mapping would not be equivariant.\n",
    "To see this, consider a vector $(1,0)^T$ which under a rotation by $\\pi$ transforms to $(-1,0)^T$.\n",
    "Acting on the former by a ReLU returns the vector itself while the latter will be projected to the origin such that all information of the feature is lost.\n",
    "\n",
    "As an alternative, we can apply nonlinearities which solely act on the *norm* of the vector but not on its orientation.\n",
    "Since the vectors rotate under the action of $C_4$, both operations commute which means that such\n",
    "[NormNonLinearities](https://quva-lab.github.io/e2cnn/api/e2cnn.nn.html#normnonlinearity)\n",
    "are equivariant for this field type.\n",
    "For instance, we can apply so called *norm-ReLUs*:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "norm_relu = nn.NormNonLinearity(feat_type_out, 'n_relu')\n",
    "\n",
    "y = norm_relu(model(x))\n",
    "\n",
    "for g in r2_act.testing_elements:\n",
    "    x_transformed = x.transform(g)\n",
    "    y_from_x_transformed = norm_relu(model(x_transformed))\n",
    "    \n",
    "    y_transformed_from_x = y.transform(g)\n",
    "    \n",
    "    assert torch.allclose(y_from_x_transformed.tensor, y_transformed_from_x.tensor, atol=1e-5), g"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Mixed Feature Fields"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Our implementation allows for different feature field types in a feature space.\n",
    "To achieve this, one simply has to pass a list of different representations when instantiating the `FieldType`.\n",
    "As an example, let's build a feature field with 2 scalar fields, 2 regular fields and 1 vector field:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [],
   "source": [
    "feat_type_out = nn.FieldType(r2_act, 2*[r2_act.trivial_repr] + 2*[r2_act.regular_repr] + 1*[r2_act.irrep(1)])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "or equivalently"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [],
   "source": [
    "scalar_fields = nn.FieldType(r2_act, 2*[r2_act.trivial_repr])\n",
    "regular_fields = nn.FieldType(r2_act, 2*[r2_act.regular_repr])\n",
    "vector_field = nn.FieldType(r2_act, 1*[r2_act.irrep(1)])\n",
    "feat_type_out = scalar_fields + regular_fields + vector_field"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "While the convolution layers can deal with mixed field types automatically, some additional care is needed when we apply non-linearities to mixed field types since these might not all allow for the same nonlinearity.\n",
    "\n",
    "In the example above, we can apply ReLUs to the scalar and regular fields but we need to apply e.g. Norm-ReLUs to the vector field.\n",
    "This can be achieved by using\n",
    "[`MultipleModule`](https://quva-lab.github.io/e2cnn/api/e2cnn.nn.html#e2cnn.nn.MultipleModule)s\n",
    "which split the input tensor into multiple branches, associating a different label to each of them.\n",
    "Each branch is then being acted on by its own module, here an `nn.ReLU` for the `'relu'` branch and an `nn.NormNonLinearity` for the `'norm'` branch.\n",
    "The output of the module is merged back together automatically."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [],
   "source": [
    "relu = nn.ReLU(scalar_fields + regular_fields)\n",
    "norm_relu = nn.NormNonLinearity(vector_field)\n",
    "\n",
    "nonlinearity = nn.MultipleModule(\n",
    "                    feat_type_out,\n",
    "                    ['relu']*len(scalar_fields+regular_fields) + ['norm']*len(vector_field),\n",
    "                    [(relu, 'relu'), (norm_relu, 'norm')]\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We once again verify that this operation is $C_4$-equivariant:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [],
   "source": [
    "model = nn.SequentialModule(\n",
    "    nn.R2Conv(feat_type_in, feat_type_hid, kernel_size=3),\n",
    "    nn.InnerBatchNorm(feat_type_hid),\n",
    "    nn.ReLU(feat_type_hid),\n",
    "    nn.R2Conv(feat_type_hid, feat_type_out, kernel_size=3),\n",
    "    nonlinearity,\n",
    ").eval()\n",
    "\n",
    "x = torch.randn(1, 1, 17, 17)\n",
    "x = nn.GeometricTensor(x, feat_type_in)\n",
    "\n",
    "y = model(x)\n",
    "\n",
    "for g in r2_act.testing_elements:\n",
    "    x_transformed = x.transform(g)\n",
    "    y_from_x_transformed = model(x_transformed)\n",
    "    \n",
    "    y_transformed_from_x = y.transform(g)\n",
    "    \n",
    "    assert torch.allclose(y_from_x_transformed.tensor, y_transformed_from_x.tensor, atol=1e-5), g"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
